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CONSTITUTIVE  RELATIONS  FOR  A  PLANAR, 
SIMPLE  SHEAR  FLOW  OF  ROUGH  DISKS 


Hayley  H.  Shen  and  Mark  A.  Hopkins 


INTRODUCTION 

The  flow  of  a  granular  material  is  an  important  transport  process  in  industry 
as  well  as  in  nature.  Examples  of  granular  flows  are  grain  flows  in  chutes  and  hop¬ 
pers,  slurry  flows  in  pipelines,  avalanches,  and  the  transport  of  sediment  and  broken 
ice  floes  in  rivers  and  oceans.  At  present,  engineers  who  work  with  such  flows 
depend  largely  upon  empirical  or  semi-empirical  information.  Recent  developments 
in  the  theories  describing  the  constitutive  relationships  for  rapidly  sheared  granular 
flows  have  revealed  many  insights  into  the  properties  of  these  flows.  These  develop¬ 
ments  may  soon  lead  to  a  sound  theoretical  basis  for  the  engineering  design  of 
granular  flow  systems  as  well  as  provide  a  better  understanding  of  many  geophys¬ 
ical  phenomena. 

However,  there  is  a  major  defect  in  the  existing  theories  that  are  capable 
of  analytically  describing  stresses  in  rapidly  sheared  granular  flows:  these  theories 
neglect  the  rotational  motion  of  particles.  From  a  laboratory  study  of  gravity 
flow  of  rough  disks  in  a  frictionless  chute  with  rough  walls,*  we  have  observed 
particle  rotation  in  the  shear  field  and  found  it  to  be  an  important  component  in 
the  interactions  between  disks.  A  sequence  of  pictures  taken  from  the  laboratory 
study  illustrates  this  phenomenon  (Fig.  1). 

Recent  computer  simulations  of  simple  shear  flows  of  rough  disks  in  a  channel 
with  rough  walls  also  showed  the  importance  of  the  particle  rotation  lu  the  rheologi¬ 
cal  behavior  of  fast  granular  flows.  Campbell  (1982)  and  Walton  (1984)  found  that 
the  rotation  of  the  disks  could  contain  energy  as  high  as  20%  of  the  total  energy 
of  the  disks'  motion. 

There  are  a  few  studies  that  have  considered  the  rotation  of  particles  (Kanatani 
1979,  Shahinpoor  and  Siah  1981).  In  these  studies,  the  rotation  is  assumed  to  be 
an  independent  kinematic  variable  unrelated  to  its  generating  mechanisms,  which 
include  the  material  properties  and  the  gradient  of  the  translational  velocity  field. 

As  a  result,  only  very  preliminary,  qualitative  expressions  for  the  constitutive  rela¬ 
tionships  are  obtained,  which  are  of  limited  engineering  interest. 

Recently,  Lun  and  Savage  (in  press)  analyzed  the  effect  of  particle  rotation 
on  the  constitutive  relations  for  a  granular  flow  of  spheres.  They  adopted  the  kinetic 
theory  approach  and  analyzed  the  complete  dynamics  of  colliding  spheres.  They 
quantified  the  rotational  motion  by  using  the  tangential  component  of  the  relative 
velocity  and  the  surface  property  of  the  colliding  spheres.  The  surface  property 
of  the  colliding  spheres  was  assumed  to  be  a  material  constant  B,  -1  s  8  £  1» 
such  that  in  the  tangential  direction  of  collision,  the  relative  velocity  before  and 
after  the  collision  has  a  ratio  equal  to  B .  The  value  B  =  1  corresponds  to  a  smooth 


*  Part  of  an  ongoing  experimental  study  at  Clarkson  University  under  National 
Science  Foundation  Grant  No.  CME-8011601. 
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Figure  1.  Time  lapse  image  of  a  gravity-driven  disk  flow  (At  =  0.04  sec). 

surface  with  no  friction,  0=0  corresponds  to  infinite  friction  and  6  =  -1  corres¬ 
ponds  to  a  complete  reversal  of  tangential  velocity.  As  Lun  and  Savage  (in  press) 
stated,  this  parameter  B ,  instead  of  being  a  constant,  should  be  a  function  of  the 
frictional  coefficient,  the  normal  impulsive  force  during  collision,  and  the  tan¬ 
gential  inelasticity  at  the  contact  surface. 

In  this  report,  we  present  an  analysis  that  quantifies  the  constitutive  rela¬ 
tions  for  a  planar,  simple  shear  flow  of  rough  disks.  This  type  of  flow  may  be  used 
to  model  the  motion  of  floating  broken  ice  in  northern  rivers  and  arctic  seas; 
therefore  the  result  is  of  practical  interest.  Also,  there  are  computer-simulated 
data  (Campbell  1982)  for  this  type  of  flow  that  can  be  used  for  direct  compari¬ 
son. 

The  rotational  motion  of  the  shearing  disks  is  computed  from  a  statistical 
average  of  a  random  collision  process  in  the  granular  material.  These  collisions 
are  inelastic  and  frictional.  It  is  assumed  that  all  colliding  disks  separate  with 
no  relative  tangential  velocity.  This  assumption  enables  us  to  directly  compare 
our  results  with  Campbell's  data.  However,  this  assumption  would  mean  an  in¬ 
finite  friction  coefficient.  Extension  of  this  analysis  to  incorporate  a  finite  fric¬ 
tion  coefficient  is  currently  underway. 

The  random  collision  process  obtains  its  energy  from  the  translational  mean 
velocity  field.  These  collisions  produce  a  random  motion  in  the  flowing  mater¬ 
ial.  This  random  motion  has  both  translational  and  rotational  modes.  Both 
modes  of  motion  are  completely  determined  by  the  imposed  mean  translational 
velocity  field  and  the  material  properties.  By  introducing  the  restitution  and 
frictional  coefficients  as  independent  parameters,  the  relative  motion  of  collid¬ 
ing  particles  in  the  tangential  direction  can  be  determined  as  a  function  of  not 
only  the  frictional  property,  but  also  the  normal  impulsive  force  at  contact.  This 
approach  incorporates  slightly  more  fundamental  physics  than  that  of  Lun  and  Savage 
(in  press).  ~ 


MODELING  PROCEDURE 


Consider  an  assembly  of  shearing  uniform  disks  as  shown  in  Figure  2,  where 
the  diameter  of  disks  and  the  size  of  the  average  gap  between  disks  are  denoted 
by  D  and  s.  The  thickness  of  the  disks  is  d,  the  density  of  the  disks  is  ps  and  the 
interstitial  fluid  is  neglected.  The  area  concentration  C  of  the  disks  is  related 
to  D  and  s,  as  shown  in  eq  1,  through  the  densest  possible  concentration  Cq  (Shen 
1982).  In  this  study,  Cq  is  assumed  to  be  0.906,  the  concentration  for  a  hexagon¬ 
al  close  packing. 


(1) 


The  disks  in  Figure  2  have  a  restitution  coefficient  e  and  friction  coefficient  p. 

In  our  study,  p  =  00  is  adopted  in  order  to  model  collisions  such  that  relative  mo¬ 
tion  at  the  contact  point  will  be  destroyed.  This  assumption  enables  us  to  make 
a  direct  comparison  with  existing  computer-simulated  data  (Campbell  1982)  and 
it  greatly  simplifies  the  analysis. 

There  are  four  kinematic  variables  for  each  individual  disk:  ui,  the  mean 
translational  velocity;  v',  the  fluctuation  of  linear  velocity;  ftk,  the  mean  angular 
velocity;  and  w',  the  fluctuation  of  angular  velocity.  Both  the  linear  and  angular 
fluctuations  are  random  variables  that  have  a  spectrum  of  magnitudes.  However, 
they  will  be  simplified  here  so  that  the  instantaneous  velocities  of  a  disk  are 


V  =  ui  +  v' 


and 


u  =  ftk  ±  a)' 


(2a) 


(2b) 


where  the  two  fluctuations  v'  and  w'  have  constant  magnitudes  and  random  direc¬ 
tions.  The  translational  velocity  v',  which  occurs  in  the  xj,  Xg  plane,  can  take  all 
possible  directions  within  the  angle  2ir ,  while  the  rotations  u>'  must  all  be  directed 
along  the  k  coordinate  axis.  Therefore 

«' =  ±w'k.  (3) 


The  magnitude  of  the  three  kinematic  variables  v',  u>'  and  ft  are  determined  by  the 
material  properties  of  the  disks  and  the  driving  mechanism  created  by  the  shear 
motion  as  defined  by  u(xg).  Without  this  shearing,  collisions  can  not  be  maintained 

in  the  flow,  hence  there  will  be 
no  established  random  motion  in 
either  the  translational  or  the 
rotational  sense. 

The  mean  angular  velocity 
ft  can  be  obtained  through  the 
following  analysis.  Consider  two 
colliding  disks  as  shown  in  Figure 
3.  Since  fluctuations  have  zero 
net  effect  on  the  averaged  rota¬ 
tion,  these  two  disks  can  be  as¬ 
sumed  to  move  with  the  mean 
velocities  as  shown  in  the  figure. 
The  relative  motion  of  the  two 


1 

1 
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Figure  3.  Colliding  disks  moving  with 
mean  velocities. 


average  sense.  Since 


points  on  the  surfaces  of  two  disks  at 
contact  is 

Uc  =  fty  +  sin1 6  +  ft  f  .  (4) 

Whenever  Uc  is  not  zero,  a  frictional 
force  acts  to  retard  the  disks’  relative 
motion  at  the  contact  point.  Since  the 
mean  linear  velocity  is  maintained  by 
the  shearing  motion,  the  retardation 
contributes  to  readjusting  their  angu¬ 
lar  velocity  only.  Eventually,  as  steady 
state  is  reached,  Uc  must  vanish  in  the 


2it 


f  ^  sin28d0  =  j 


the  average  of  eq  4  yields 

Uc  =  ft  D  +  ^Dp-  =  0 
c  2  dX2 


(5) 


which  implies 


ft  =  - 


J_  du 
2  dx2 


(6) 


This  same  result  has  been  shown  by  Campbell  (1982)  in  the  computer-simulated 
shear  flow  of  disks  and  has  been  derived  by  Lun  and  Savage  (in  press)  for  a  sim¬ 
ple  shear  flow,  assuming  that  random  motions  have  a  Maxwellian  distribution. 

As  the  derivation  of  ft  shows,  when  ft  =  -|(du/dx),  only  in  the  average  sense  is 
there  no  relative  motion  between  two  colliding  disks.  Both  the  randomness  of 
contact  location  and  relative  impact  velocity  ascribable  to  disk  fluctuation  can 
cause  the  angular  velocity  to  deviate  from  the  mean  value  ft.  Denoting  the  angu¬ 
lar  velocity  fluctuation  by  w\  this  fluctuation  is  a  consequence  of  random  contact 
location  and  the  fluctuation  of  the  linear  velocity  v'. 

The  angular  velocity  fluctuation  w'  is  determined  in  the  following  way. 

The  assumption  of  steady,  simple  shear  flow  implies  that  the  energy  contained 
in  the  rotational  mode  must  be  constant.  Let 


E r  =  to'2  .  (7) 

Er  being  constant  requires  that,  on  average,  the  angular  fluctuations  before  and 
after  collision  must  be  the  same.  Consider  two  colliding  disks  A  and  B  as  shown 
in  Figure  4.  The  angular  velocities  of  these  two  disks  are  randomly  chosen  as 


toft  =  ft  ±  to' 

(8a) 

0)0  =  ft  ±  to'. 

(8b) 

Let  ui*  denote  the  post-collision  velocity  of  disk  B  from  one  of  the  four  random 
combinations  of  w ft  and  to  0.  The  condition  that  Er  =  w'2  is  constant  can  be 


Figure  4.  Colliding  disks  with  ve¬ 
locity  fluctuations. 

described  as 


Figure  5.  Energy  balance  in  a  con¬ 
trol  area. 


(0'J  = 


/  J-y  da  dQ 

f  f  da  dQ 


where  9  and  a  are  defined  in  Figure  4  and  the  integration  is  over  all  possible  val¬ 
ues  of  these  two  variables.  This  integration  will  be  carried  out  later  when  ranges 
for  0  and  a  are  determined.  The  result  of  this  integration  determines  o>’. 

The  two  kinematic  variables,  v'  and  w',  determine  the  energetics  of  this 
shearing  assembly  of  disks.  Consider  a  control  area  as  shown  In  Figure  5,  where 
W  represents  the  rate  of  work  done  by  the  shear  forces  and  $  represents  the  rate 
of  dissipation  from  inelastic  and  frictional  collisions.  In  a  homogeneous  flow,  as 
the  shearing  reaches  a  steady  state,  IV  must  be  equal  to  $.  In  mathematical  form, 
this  can  be  written  as  IV  =  $  or 


du  . ,  v'  .  =, 

T2iajr=N  t4E 


where  x  2\  is  the  shear  force,  N  is  the  number  of_disks  in  a  unit  control  area,  s  is 
the  mean  gap  size  between  adjacent  disks,  and  AE  is  the  average  colllsional  energy 
loss. 

The  stresses  are  modeled  as  a  result  of  coliisional  momentum  transfer 
across  control  surfaces  (Bagnold  1954) 


Tzi  =  P-yAM  l 


(11a) 


T  22  =  P  AMg 


(lib) 


where  p  is  the  number  of  disks  per  unit  length  along  the  boundary  of  the  control 
area  in  Figure  5,  and  AM  j  and  AMg  are  the  average  linear  momentum  transfer  in 
the  xj  and  xg  directions  because  of  collisions. 

As  6hown  by  Shen  (1982), 


p  =  ttD 


l  V-CU*  - 


and 


— r —  (13) 

C* 

where  C  is  the  area  concentration  and  Cq  is  the  densest  area  concentration  that 
corresponds  to  s  -  0.  Substituting  eq  11  into  eq  10  and  using  p/D  =  N  (Shen  1982), 
we  obtain 


s  =  D 


dti 
dx 2 


D AMi  =  HE . 


(14) 


The  quantities  AMj  and  HE  are  functions  of  the  material  properties  of  the  disks 
and  the  four  kinematic  variables  du/dxj,  ft,  v'  and  a/.  With  u>'  determined  in  eq  9, 
eq  14  is  then  used  to  solve  for  v'.  After  v'  and  A  Mi  are  determined,  the  stresses 
can  then  be  specified  using  eq  1  la  and  b. 

_  In  the  next  section,  details  are  provided  regarding  the  solutions  for  to',  v\ 

A  Mi  and  A  M2. 


DERIVATION  OF  STRESSES 

In  this  section,  the  four  averaged  quantities— angular  velocity  fluctuation 
(to1),  translational  momentum  transfer  in  the  Xi  and  X2  directions  (A  Mi  and  A  M2) 
and  energy  loss  during  collisions  (AE)— are  first  determined  as  functions  of  the 
material  properties  and  the  strain  rate.  Equations  11a  and  b  are  then  used  to 
compute  the  shear  and  normal  stresses. 

As  shown  in  Figure  6,  disks  A  and  B  are  colliding  at  the  location  specified 
by  0 .  The  coordinate  system  is  fixed  on  disk  B  and  moves  with  the  mean  velocity 
of  the  disk  just  prior  to  collision.  The  relative  fluctuation  velocity  between  A 
and  B,  / 2v',  Is  randomly  distributed  between  the  angles  a  1  and  a  z.  The  four  quan¬ 
tities  w',  AMi,  AM2  and  HE  are  obtained  by  Integrating  over  0  and  a. 

The  values  of  a  1  and  a  2  are  derived  from  a  kinematic  constraint  on  v',  such 
that  the  pre-collision  velocity  of  A,  V^g,  must  be  within  ±  it/ 2  of  N  in  order  for  a 
collision  to  occur  (see  Fig.  6).  From  the  law  of  sines,  it  can  be  shown  that  (Shen 

1982) 

a 2  =  -cx  1  =  ^  -  sin_1(B  sin20)  (15) 


where 

R  =  d£t~/  2/2v\  (16) 

uXg 

Using  the  approximation  sin-1x  ~  x 
yields 

a  2  =  -  a  1  -  — 2  -R  sin  20.  (17) 

The  unit  vectors  N  and  P  as  shown  in 
Figure  6  are  the  normal  and  tangential 
vectors  at  a  g.ven  contact  point  speci¬ 
fied  by  0 .  The  components  of  N  and  P  are 
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Figure  6.  Definition  sketch  of  ran¬ 
domly  colliding  pairs  of  disks. 


N  =  -  cos6  i  -  sine  j 


(18a) 


P  =  sine  i  -  cos6  j. 


(18b) 


The  velocity  components  of  disk  A  prior  to  the  collision  are  shown  in  Appen¬ 
dix  A  (eq  A4a  and  A4b)  as 


VA  =-D  .  sine  cose  +  /2  v'  cosa 


%  =  D^7sin*e  +  v5v'  sina. 


Average  momentum  transfer  in  the  x\  direction— A  Mj 

A  Mi  is  the  average  transfer  of  momentum  in  the  Xi  direction  attributable 
to  collisions.  By  use  of  the  reference  frame  that  moves  with  the  mean  motion 
of  disk  B  just  prior  to  collision,  the  transfer  of  momentum  from  disk  A  to  disk  B 
becomes  equal  to  the  post-collision  momentum  of  disk  B  in  the  Xj  direction. 
Therefore 


A  Mi  =  m  V*-  i 


where  V*  is  the  post-collision  velocity  of  particle  B.  The  average  momentum 
transfer  AAfj  resulting  from  all  possible  collisions  can  be  expressed  as 


"D-df'i”2  AM, da 

AAfi  =  ps  - - - 

v  a  7 

4^o4,d»  de 


where  a  i  and  a  2  are  defined  in  eq  15.  The  velocity  V*can  be  described  in  terms 
of  the  components  in  the  N  and  P  directions: 


V*  =  V*  N  +  V*  P. 
B  Bf|  Bp 


The  post-collision  velocities  of  disks  A  and  B  depend  on  the  pre-collision 
velocities  and  the  material  properties  characterized  by  the  restitution  coefficient 
e  and  the  friction  coefficient  p.  In  this  analysis,  the  post-collision  tangential 
velocity  between  A  and  B  is  assumed  to  be  zero.  This  is  the  same  condition  used 
in  Campbell's  computer-simulated  shear  flow  of  disks  (1982)  and  implies  that  the 
friction  coefficient  p  is  infinite. 

Letting  p  equal  infinity,  the  velocity  components  of  disk  B  after  the  colli¬ 
sion  are  derived  in  Appendix  A  (eq  A5c  and  A5d)  as 


.  (llgj  (_D  du. 

2  '  udx2 


sin0  cos0  +  /2v'cosa) 


V*  =  jr  {  D  ~t~  sin*0  +  / 2  v'  sina  +  -y  (a>^  +  u)g)} . 


m 


The  angular  velocities  oi*  and  tog  contain  both  mean  and  fluctuating  com¬ 
ponents,  as  defined  in  eq  2b.  However,  since  the  average  of  the  fluctuating 
components  is  zero,  aM  does  not  contain  the  fluctuation  component  of  angular 
velocities.  Hence  only  the  average  angular  velocities  are  necessary  to  describe 
Vjjfp.  From  eq  6  it  is  seen  that 

n  1  du 

“A=  “B=  0  =  -  2  d3F7 

therefore,  V£p  can  be  written  as 

%  cos20+v'2v'sina).  (22) 

After  substituting  eq  18a  and  b,  21  and  A5c  and  d  into  eq  20  and  integrat¬ 
ing,  we  obtain  the  average  momentum  transfer  in  the  Xj  direction  as 


^  =P s  d(D-|p)(0.106e  +  0.124). 


(23) 


Average  momentum  transfer  In  the  direction^ M2 

A  M2  is  the  average  transfer  of  momentum  in  the  X2  direction  ascribable  to 
disk  collisions.  Again  our  reference  frame  is  fixed  on  disk  B  prior  to  a  collision. 
In  this  reference  frame,  the  momentum  transferred  is  equal  to  the  post-collision 
momentum  of  disk  B  in  the  X2  direction.  Hence  the  average  value  of  the  momen¬ 
tum  transfer  in  the  X2  direction  is 


am2  = 


ttDM  Jj^^-jda  d6 
~~~A  /"J^dode 

Jo  Jai 


(24) 


with  aj  and  012  defined  in  eq  15. 

After  substituting  eq  1 8a  and  b,  21  and  A5c  and  d  into  eq  24  and  integrat¬ 
ing,  we  find  that  the  average  momentum  transfer  in  the  x2  direction  becomes 


A M2=  -Ps1^  {0.012(D^)2  +  0.287(l+e)v'2}.  (25) 


Rotational  fluctuation-id ' 

The  instantaneous  angular  velocity  oi  has  both  a  mean  and  a  fluctuating 
component  as  defined  earlier  in  eq  2b.  Since  all  rotations  are  in  the  k  direction, 
w  may  be  written  as 

0)  =  ft  ±  oi'.  (26) 


Consider  two  colliding  disks  A  and  B  as  shown  in  Figure  6.  Let  oi*  denote 
the  post-collision  angular  velocity.  As  shown  in  Appendix  A  (eq  A2e  and  f  with 

8  =  i), 


WA=  w8  = 


(27) 


where  u>a  and  013  are  given  in  eq  26.  Since  there  are  four  arbitrary  combina¬ 
tions  of  oiA  and  013,  the  post-collision  angular  fluctuations  are 


i\ 


1  VAd 

o>*  -  fi  =  j  { ±2u,'  ±  w'-2(n+  )}  (28) 

where  the  four  combinations  of  the  ±  signs  are  to  represent  the  four  combina¬ 
tions  of  u  a  and  uiq. 

From  Appendix  A  (eq  A4b) 

j  VAP  =  D  Tfitj  sine +/2  v'  sina . 

The  average  angular  fluctuations  as  derived  in  eq  9  can  now  be  obtained  after 
substituting  eq  28  and  A4b  into  eq  9.  The  integration  limits  are  again  0  s  0  ^  it 
and  oil  =  “  -  oi2»  w^h  a  i ,  a2  given  in  eq  17.  This  integration  then  gives 

«*  -  0-125  (-^  ) 2  *  1.024  (£ )  \  (29) 

This  derivation  for  to1  is  rather  formal.  A  more  intuitive  way  of  deriving 
u'  is  given  in  Appendix  B,  where  a  6emi-Monto-Carlo  method  is  used  to  generate 
<*)'.  The  result,  however,  is  identical  to  eq  29. 

Energy  loss-AE 

AE  is  the  average  energy  loss  in  an  inelastic  and  frictional  collision.  As 
discussed  earlier  (eq  7),  the  average  rotational  energy  does  not  change  in  a  steady 
state.  Hence  A  E  may  be  modeled  as  the  result  of  translational  energy  loss  alone. 
For  a  given  collision,  the  translational  energy  loss  in  a  pair  of  disks  A  and  B  is 

AET  =  \  m(VA-  VA  -  V**  V*  +  VB*  VB  -  V*  •  V*)  (30) 

where  m  =  ps(ir  D  2d/4)  is  the  disk's  mass.  ___ 

The  average  translational  energy  loss  A  Ej  is 

^.-i-  C ( < * v M  (31) 

at** 

since  YB  is  zero  in  the  chosen  reference  frame.  As  derived  In  Appendix  A  (eq 
A5a  and  b). 


V*  =ii±Li  Sjn0  COS0  +  /2v'cosa) 

Afj  2  ax2  1 

v *  =  )(5  sln20+|  )  +  5/2  v' sina 


Substituting  eq  A4a  and  b  and  A 5a,  b,  c  and  d  into  eq  31  yields  the  average  colli- 
sional  energy  loss 

AE  =  (0. 0464-0. 0261e^)  +  v'2  (0. 386-0. 244e2)  .  (32) 
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In  a  steady,  homogeneous  shear  flow,  the  work  done  by  shear  forces  must 
equal  the  energy  dissipation  from  inelastic  and  frictional  effects  of  the  colli¬ 
sions  between  disks.  This  balance  law  has  been  shown  to  reduce  to  (eq  14) 

=  A E. 

Substituting  eq  23  and  32  into  eq  14,  we  obtain  the  fluctuation  of  linear 
velocity 


.  ..  du  a  ,  0.0261c2  +0.106e  +0.0776  V2 

V  - 0.386-0.244? - >  ' 


(33) 


Substitution  of  eq  33  into  eq  29  allows  the  flucutation  of  angular  velocity  to  be¬ 
come 


,  _  _du_  ,  -0.0038E2  +  0.109e +  0.128  ^2 
“  dx2  1  0.372-0.244^  J 


(34) 


At  this  point  all  the  kinematic  variables  have  been  found  in  terms  of  the 
mean  velocity  gradient. 


SHEAR  AND  NORMAL  STRESSES 

The  shear  and  normal  stresses  are  obtained  by  substituting  eq  23,  25,  and 
33  into  eq  1  la  and  b.  This  yields 

'a  ■  ^(°)'D«»2(^Aoao6£+o.u4)(«5«^)1/2  ,35*) 


and 


du  0.0075e3  +  0.035e2  +  0.053e  +  0.027  /,ruX 

PsC(T)(DM)2(-5— )  - 0.386-0.244^ -  (35b) 


dx2 

Removing  the  dimensions  from  the  stresses  by  dividing  eq  35a  and  b  by 
PsCq D^du/dx^V  results  in  the  following  nondimensional  stresses, 


(0.106e+0.124)  r  0.0261c2  +  0.106e  +  0.0776^ 1/2 
T21  ~  [(Cq/C)1^-!]  ^  0.386-0. 244^  > 


(36a) 


and 


t22  =  Rc07c)^-i] 


-1 


0.0075e3  +  0.035c2  +  0.053e  +  0.027 
0.386-0.244C2 


(36b) 


COMPARISON  WITH  COMPUTER-SIMULATED  STRESSES 

To  date,  there  have  been  no  data  obtained  from  physical  models  of  shear¬ 
ing  disks  that  can  be  used  to  verify  the  theory  developed  above.  The  only  ex¬ 
perimental  data  currently  available  are  from  computer  simulations  by  Campbell 
(1982).  The  comparison  between  Campbell's  data,  eq  40a  and  b  and  the  shear 


c 


c 


Figure  7.  Comparison  of  theoretical  prediction  and  computer 
simulation. 


stress  obtained  by  Shen  (1982),  which  did  not  consider  disk  rotation,  is  presented 
in  Figure  7. 

We  observe  that  by  including  disk  rotations,  the  shear  stress  is  reduced 
for  the  case  when  y  =  ».  The  comparison  between  the  computer  simulation 
data  and  eq  36a,  b  is,  in  general,  satisfactory  except  when  the  disk  concentra¬ 
tion  is  lower  than  0.3.  We  believe  that  at  lower  concentrations,  the  momentum 
transfer  that  is  related  to  mass  transfer  increases.  However,  in  the  present 
analysis,  mass  transfer  is  not  considered.  Since  for  high  concentrations,  mass 
transfer  is  negligible,  the  correlation  is  closer  in  this  regime.  As  the  concentra¬ 
tion  decreases,  we  therefore  expect  that  collisional  momentum  transfer  alone 
would  produce  stresses  that  are  lower  than  measured  values. 


DISCUSSION 

We  make  two  major  observations  on  the  result  of  this  theoretical  develop¬ 
ment.  The  first  is  that  the  assumption  of  equipartition  of  fluctuation  energy 
does  not  apply  in  a  dissipative  granular  flow.  The  second  is  that  the  distribution 
function  of  the  random  fluctuating  states  does  not  have  a  significant  effect  on 
stresses. 

Since  the  modeling  procedure  for  determining  the  stresses  is  very  similar 
to  the  molecular  dynamic  approach  used  with  dense  gasses,  and  since  the  equi¬ 
partition  of  energy  in  molecular  random  motion  is  well-known  for  gasses  that 
are  in  equilibrium,  it  has  been  suggested  that  a  form  of  equipartition  between 
the  rotation  and  translation  motion  of  each  individual  disk  may  exist  (Kanatani 
1979).  It  is  therefore  interesting  to  compute  the  ratio  between  fluctuation 
energy  in  the  rotational  motion  and  that  in  the  translational  motion,  which 
yields 
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WUJULW4HJH  JVM.” U1  1 JU  JI 1  JlALfclllU  UV 


h  -  tfaj'2.  (37) 

rt  "  fmv72 

where  I  is  the  moment  of  inertia  for  disks.  Substituting  eq  29  and  33-and  using 
I/m  =  D2/ 8  for  a  disk-lnto  eq  37  gives  us 


0.06-0. 04c2 
0.3c2  +  l.le  +  0.8 


■)  +  0.128 


(38) 


which  is  far  from  being  equal  to  unity  and  depends  explicitly  on  the  material 
property  e .  The  reason  that  the  friction  property  does  not  appear  in  eq  38  is 
that  the  condition  of  y  approaching  infinity  is  assumed  in  this  analysis.  We  ex¬ 
pect  that  when  a  finite  friction  coefficient  is  considered,  the  partition  ratio 
would  also  be  a  function  of  the  friction  coefficient. 

In  this  theoretical  development,  it  is  recognized  that  both  the  linear  and 
angular  velocities  of  disks  are  random  variables.  However,  in  order  to  make  the 
mathematical  manipulation  manageable,  serious  simplifications  have  been  made. 
We  assumed  that  the  fluctuation  components  v'  and  y’  are  randomly  directed, 
but  constant  in  magnitude.  In  order  to  make  a  constant  v’  and  w’  contain  the 
same  energy  as  a  distributed  v'  and  y',  the  constant  v'  and  y'  must  be  the  stan¬ 
dard  deviations  of  the  real  distributions.  It  has  been  suggested  that  the  calcu¬ 
lated  stress  level  could  be  lowered  significantly  by  ignoring  the  real  distributions 
(which  could  be  Maxwellian  as  suggested  by  Jenkins  and  Savage  [1983]  and  Lun 
et  al.  [1984])  and  by  assuming  a  constant  magnitude  for  the  fluctuating  compo¬ 
nents.  The  stresses  obtained  by  Shen  and  Ackerman  (1982),  for  example,  when 
analyzing  the  simple  shear  flow  of  uniform  spheres  were  about  three-fold  lower 
than  those  obtained  by  Lun  et  al.  (1984). 

However,  the  result  of  the  present  analysis  suggests  that  this  is  not  the 
case,  at  least  for  a  flow  of  disks. 

This  conclusion  is  based  upon  the  close  agreement  between  the  present 
theoretical  results  and  the  computer-simulated  data  of  Campbell  (1982),  since 
in  the  computer  simulation,  the  distribution  of  fluctuation  motion  is  inherently 
included.  If  the  velocity  distribution  can  make  a  three-fold  magnification  of 
stresses,  it  would  also  be  revealed  in  Campbell's  data.  The  fact  that  the  present 
theory  and  Campbell's  data  agree  so  well  supports  the  validity  of  this  simple 
model. 


CONCLUSIONS 

Constitutive  relations  are  developed  for  a  granular  continuum  undergoing 
a  simple  shear  flow.  The  granular  continuum  flow  consists  of  uniform  disks  mov¬ 
ing  in  a  plane.  The  result  of  the  analysis  may  be  used  to  model  the  motion  of 
broken  ice  fields  on  rivers  and  oceans. 

The  major  contribution  of  this  analysis  is  that  it  includes  the  disk  rotations. 
The  analysis  was  able  to  relate  this  rotation  entirely  to  its  generating  mechanisms, 
namely,  the  linear  velocity  gradient.  Previous  studies  of  this  problem  either  ig¬ 
nored  the  rotation  completely  (e.g.  Shen  and  Ackermann  1982,  Jenkins  and  Sav¬ 
age  1983),  or  treated  rotation  as  an  independent  mode  of  motion,  unrelated  to 
the  linear  velocity  gradient  (e.g.,  Kanatani  1979,  Shahinpoor  and  Siah  1981). 

The  stresses  are  formulated  as  the  rate  of  momentum  transfer  attributa¬ 
ble  to  binary  collisions.  These  collisions  are  inelastic  and  frictional.  The  re¬ 
sults  show  explicitly  the  effect  of  material  properties  on  the  stress  level.  How¬ 
ever,  the  constitutive  laws  remain  in  the  second  power  form 


T 


_duj)2 
dx  2 


which  makes  granular  flow  behave  differently  from  viscous  fluid  flows. 

The  significant  findings  of  this  investigation  are: 

1.  The  stress  level  is  decreased  from  previous  analysis  (Shen  and 
Ackermann  1982)  when  disk  rotating  is  included. 

2.  There  is  no  equipartition  of  rotation  energy  and  translational 
energy. 

3.  The  distribution  of  random  fluctuation  does  not  have  significant 
effect  on  stresses. 

Good  agreement  between  the  theoretical  results  and  computer-simulated 
data  is  obtained.  This  demonstrates  the  need  for  continued  development  of  the 
simple  model  used  in  this  analysis.  Future  work  should  extend  the  present  model 
to  analyze  mixtures  of  different  sizes  and  irregular  shapes.  In  order  to  analyze 
a  general  flow,  the  analysis  should  also  be  extended  to  include  a  general  deforma¬ 
tion  field  described  by  a  tensor  with  all  four  components  Suj/Sxi,  3ui/3x2, 
3u2/3ui,  3u2/3u2.  The  above-mentioned  studies  are  currently  underway. 
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APPENDIX  A:  POST-COLLISION  VELOCITIES 


I  Consider  two  disks  A  and  B  colliding  with  each  other  as  shown  in  Figure  6. 

i  Equations  governing  the  collision  between  disks  A  and  B  were  developed  using 

I  the  laws  of  conservation  of  linear  and  angular  momentum,  the  energy  loss  be- 

[  cause  of  incomplete  restitution,  as  specified  by  the  coefficient  of  restitution, 

i  and  a  frictional  closure  condition. 


(Ala) 

%  +  %  =  VVBp 

(Alb) 

vAp  % 

8wa  r  -  r 

(Ale) 

vAp  ,  % 

B<*>B+_r~  =  8wb+-^- 

(Aid) 

<%  -  /  <%-  V  =  ~  e 

(Ale) 

(A  If) 

The  above  six  equations  are  used  to  solve  for  the  six  post-collision  varia¬ 
bles,  V/JL  VgJj,  Vtf  ,  w*  and  t»*.  The  frictional  closure  condition  used  in  this 
report  is  that  the  relative  tangential  velocity  at  the  point  of  contact  after  the 
collision  is  zero.  This  implies  that  friction  is  always  great  enough  to  destroy 
the  relative  tangential  velocity,  or  that  friction  coefficient  is  infinite.  The  six 
post-collision  variables  are  obtained  from  eq  Ala-f  as 


%  -  '-r1  % 
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(A2a) 


<B  +  2)  VA  +  B  VBp  -DB (u^v+Oj  b)  /  2 
=  2TU6l 
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(1+e) 
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(A  2b) 

(A2c) 
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%  =  2(ITB) 


(A2d) 


_  (1+2B)  oja-wr-2  <typ-VBp)  /  D 
2(1+B) 

-wa+(1+2B)wB-20 Ap-VBp)/D 
“b  = - 2(1+3)" 


(A2e) 

(A2f) 


B  is  the  ratio  of  the  square  of  the  radius  of  gyration  to  the  square  of  the 
disk  radius.  For  a  uniform,  circular  disk,  B  =  j. 

With  the  reference  frame  located  at  the  center  of  disk  B  prior  to  a  colli¬ 
sion  as  shown  in  Figure  6,  the  pre-collision  velocity  of  disk  B  is  VB  =  0.  There¬ 
fore,  eq  A2a-f  can  be  simplified  to 
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(A3b) 
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(A3e) 

>  r  *»  2V‘-  i 

=  3  [-wa+2<»)b - YT  )' 


From  Shen  (1982) 

V/i,,  =  -  D  Sin0  COS0  +/ 2  v'  cosa 
"N  ax  2 

VAp  =  D-^  sin20+/2  v'  sinot  . 

Substituting  these  into  eq  A3a-f,  we  obtain 

=  ^2  ~  [~D-jl~slnQ  cos0+/2v'  cosa) 

V*^  =  ^  (5D-~sin20-5/2v'  sina  --|(uA+u)B)) 


V*  _  (l±Ll  (-d~-  sin0  cos0  +  /2  v'  cos  a  ) 

2  '  ax  2 

^Bp  =  ^  (D-^jsin20  +  /2  vf  sin  a+ (ua+wb)) 

a»A  =  wB*  =  “  j  (|  +  2  sin20)  +  2/2  (v'/D)  sina  . 

In  the  process  of  averaging  the  following  series  expansions  were  used: 

x3  x5  x7 

sinx  =  x-  3r+5j--7r+... 

r2  r6 

cosx  =  1  -  2f  +  4J  -  gT  +  •  •  • 

where  x  =  a2  =  ir/2  -  R  sin  20  and  R  =  D(du/dx2)/ 2/2  v'.  The  series  were  trun¬ 
cated  after  the  Ji2  term. 


du 


(A3f) 

(A4a) 

(A4b) 

(A5a) 

(A5b) 

(A5c) 

(A5d) 

(A5e) 
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APPENDIX  B:  AN  ALTERNATIVE  METHOD  OF  FINDING  o»’ 


We  have  defined  the  post-collision  angular  velocity  of  a  disk  to  be 


to*  =  12  ±  to'.  (Bl) 

Using  eq  A3e,  we  can  define  a  "second  generation"  post-collision  angular 
velocity  to  be 


2VAb 


0) 


**  1  *  *  "P 

A  =  “  «B - )  • 


(B2) 


Substituting  eq  Bl  into  eq  B2,  and  subtracting  the  mean  component,  ft, 
yields  a  "second  generation"  fluctuating  angular  velocity,  which  we  then  square 


(B3) 


If  we  begin  with  oi’  =  oig  =  0  (first  generation),  we  find  to  be  a  function 
of  (ft+V^p/D),  hence  each  succeeding  generation  of  to'  is  also  a  function  of 
(ft+  VAp/D),  which  can  thus  be  moved  outside  the  brackets.  This  quantity  is  a 
function  of  the  point  of  collision  denoted  by  Gand  the  direction  of  v'  denoted  by 
a  and  is  averaged  by  integration  over  the  range  0  and  a. 

We  iterate  eq  B3  over  all  previous  values  of  to'  using  the  following  expres¬ 
sion 
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%  ,2 
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(B4) 


where  w.'  =  to'  =  0  for  f  =  /  =  0. 

As  n  grows  large  (~  50)  the  bracketed  expression  approaches  unity.  Evaluat 
ing  the  integral  expression  yields 


„.2=°.,25(-^)2* 


(B5) 


which  is  iuentical  to  the  result  obtained  by  the  method  outlined  previously  in 


A  facsimile  catalog  card  in  Library  of  Congress  MARC 
format  is  reproduced  below. 
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